• Fazer afirmações sobre as características de uma população, com base em informações dadas por amostras.

• População é qualquer conjunto de elementos que você está observando

• Amostra é qualquer subconjunto desta população.

Ajustando um modelo

• Estimar parâmetros da população baseados na amostra, modelo teórico de distribuição de probabilidades

• Variável Aleatória de distribuição de probabilidades Normal com média (μ) e desvio padrão (σ).

## [1] 22.72789

## [1] 32.3065

## [1] 5.683881

A variável aleatória X é uma variável que tem um valor único (determinado aleatoriamente) para cada resultado de um experimento (lembre-se esse termo veio da teoria de probabilidade).

Exemplos:

1. Altura das árvores de uma floresta (!!)

2. Número de presas capturadas por um predador em um determinado dia;

3. Comprimento de um peixe adulto selecionado aleatoriamente.

As variáveis aleatórias podem ser discretas ou contínuas.

Função que define as probabilidades P(X) para cada valor possível da variável X.

Distribuição de Probabilidades.

• Discreta, valores possíveis 0 e 1. p - probabilidade de sucesso
  
• Jogar moedas: cara - 0, coroa - 1.
  
• Presença/ausência de uma espécie em locais amostrados

• Discreta: número de sucessos em n tentativas . p - prob sucesso.

• Número sementes predadas em cada experimento contendo 20 semente.

• Número de animais infectados em relação ao número de capturados em cada local

• Predação de girinos por odonatas, p = 0.30.
  Total girinos: 6 num lago.
• Qual probabilidades de que 0, 1, 2, 3, 5 ou todos sejam predados?

• Discreta: Probabilidade de uma série de eventos ocorrem em um período fixo de tempo, área, volume, quadrante, etc.

• Dados de contagem

• Abundância de espécies em um local

• Número de indivíduos capturados por tempo

• Número de visitas de borboletas em uma planta em 15 min: 10
• Probabilidade de 5 visitas dpois(5,10)

• Qual é a probabilidade de que um peixe capturado aleatoriamente tenha 20,15 cm ou mais? Média da população é 17,1 cm, desvio padrão 1,21 cm?

## [1] 0.005856729

"Toda soma de variáveis aleatórias independentes de média finita e variância limitada é aproximadamente Normal, desde que o número de termos da soma seja suficientemente grande."

Independentemente do tipo de distribuição da população, na medida em que o tamanho da amostra aumenta, a distribuição das médias amostrais tende a uma distribuição Normal.

• Dados de duas espécies

levels(arv$sp)

## [1] "sp1" "sp2"

Será que a altura média de uma espécie nessa floresta é diferente da altura média da outra espécie?

Ou seja, será que as espécies possuem uma mesma distribuição de probabilidades e as diferenças encontradas são devido ao acaso ou será que cada espécie é uma variável aleatória com médias diferentes?

mean(arv$alt[arv$sp=="sp2"]) - mean(arv$alt[arv$sp=="sp1"])

## [1] 5.109466

Posso afirmar que as médias são diferentes?

H0 - não há diferença entre as alturas das duas espécies. A média da altura da espécie 1 é igual à média da altura da espécie 2 e a diferença observada se deve ao acaso.

H1 - há diferença entre as alturas das duas espécies. A média da altura da espécie 1 é diferente da média da altura da espécie 2.

AFIRMAÇÃO: não há diferença entre as alturas das espécies,
VERDADE: há diferenças.

AFIRMAÇÃO: há diferença entre as alturas das espécies,
VERDADE: não há diferenças.

Definição formal dos erros:

ERRO DO TIPO I (α): rejeitar a hipótese nula dada que ela é verdadeira.

ERRO DO TIPO II (B): aceitar a hipótese nula dado que ela é falsa.

Contrastar a probabilidade de que a diferença que você encontrou nas médias das alturas das duas espécies pode se considerada devido ao acaso ou não.

Qual é a probabilidade de que as duas amostras pertencem à mesma população de medidas?

Podemos ou não rejeitar a hipótese nula (de que a diferença é ao acaso), baseado em alguma probabilidade de estarmos cometendo algum erro (α).

"O número de observações independentes menos o número de parâmetros estimados."

Quantas observações independentes foram usadas para calcular a média e a variância dos dados sob a hipótese nula?

Nosso exemplo: 100 árvores - 2 parâmetros = 98 graus de liberdade

Não entraremos em detalhes de como se calcula a estatística t e como é a sua distribuição de probabilidades (recomendo ver na wikipedia). O importante a saber aqui é que tendo as duas amostras modeladas como uma distribuição normal, calculamos a estatística t e com um certo valor de graus de liberdade, nós podemos recorrer às tabelas da distribuição t para ver qual a probabilidade de se obter aquele valor que obtivemos dada que a hipótese nula é verdadeira. Se esta probabilidade for muito pequena, a chance de estarmos caindo no erro do tipo 1 é pequena.

t.test(arv$alt~arv$sp, var.equal=T)

## 
##  Two Sample t-test
## 
## data:  arv$alt by arv$sp
## t = -5.0125, df = 98, p-value = 2.387e-06
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -7.132311 -3.086621
## sample estimates:
## mean in group sp1 mean in group sp2 
##          20.17316          25.28262

#o argumento var.equa=T diz que as variâncias das duas espécies são as mesmas
# podemos fazer o teste com a premissa de que as variâncias não sejam iguais,
# usando var.equal=F, isso vai "consumir" graus de liberdade.

Vamos então incluir no nosso exemplo inicial uma terceira espécie amostrada (e vamos diminuir o número de amostras por espécie para ficar mais fácil de fazer os cálculos passo a passo):

alt.sp3 <- c(36.90076761,33.23131727,32.82767311,24.93410577,31.87626329,
         28.76824248,26.6436144,31.41238398,27.65383929,31.70425719,
         33.60752445,28.30263354,29.32210674,23.4790054,22.20793046,
         18.28797293,23.39044736,28.75820917,29.94344703,34.47430235,
         32.09820862,26.99845592,27.64858951,31.90557306,33.63236368,
         25.35966179,35.00885172,32.89125598,18.98802229,30.8516906,
         22.04191398,35.37232889,35.67150408,34.30257037,24.43038475,
         20.6274663,31.15717916,23.27023231,36.20568545,34.33068448,
         28.40769054,31.00578874,29.50989158,31.69915991,37.25975797,
         30.1279997,38.623371,35.36472829,29.88304259,31.00145491)

sp3 <- data.frame(alt=alt.sp3, sp=rep("sp3",50))
arv2 <- rbind(arv,sp3)

amostras <- seq(1,150,by=5)

arv3 <- arv2[ amostras, ]

Neste exemplo lidaremos com um modelo simples de anova chamado one way ANOVA, porque lida apenas com uma variável independente categórica (fator).

A ANOVA testará a hipótese nula de que as médias das alturas das 3 espécies não diferem.

Faremos agora a ANOVA "na unha", para entendermos melhor cada passo da análise e a interpretação de seus resultados para depois usarmos a função do R que faz esta análise diretamente.

ANOVA na unha

Vamos calcular os valores da tabela de ANOVA. Começando com os desvios quadráticos, ou seja, quanto os dados desviam da média (idéia parecida com a variância). O ponto importante é que essa variação é aditiva e portanto, pode se decomposta. A variação total é decomposta em variação relacionada ao tratamento (entre grupos), no nosso caso às espécies, e uma variação interna dos grupos (chamada de erro). A estatística F é a razão entre essas variaões após dividir cada uma delas pelos seus respectivos graus de liberdade. Complicou? Vamos fazer os cálculos e ver os gráficos para ver se entendemos melhor:

A tabela de ANOVA que vamos construir é essa:

Primeiro vamos mudar a forma do nosso data frame para facilitar os comandos

arv4 <- data.frame(sp1=arv3$alt[arv3$sp=="sp1"],
                   sp2=arv3$alt[arv3$sp=="sp2"],
                   sp3=arv3$alt[arv3$sp=="sp3"])
arv4

##          sp1      sp2      sp3
## 1  21.282000 19.14508 36.90077
## 2  19.875790 18.66252 28.76824
## 3  23.562570 25.10581 33.60752
## 4  19.313744 18.10266 18.28797
## 5  26.308616 21.48084 32.09821
## 6   8.301795 28.08559 25.35966
## 7  27.192502 20.47805 22.04191
## 8  15.941866 29.22927 20.62747
## 9  20.596624 19.32911 28.40769
## 10 22.048945 18.21446 30.12800

boxplot(arv4)

Vamos calcular a média geral, as médias para cada espécie, e as diferenças entre a média geral e a média para cada espécie, para chegar ao valor da soma dos quadrados total:

media.geral <- mean(arv3$alt)
media.geral

## [1] 23.28284

medias.sps <- apply(arv4, 2, mean)
medias.sps

##      sp1      sp2      sp3 
## 20.44245 21.78334 27.62274

dif.geral <- arv4 - media.geral
ss.especies <- dif.geral^2
ss.especies

##                sp1       sp2        sp3
##  [1,]   4.00337202 17.121094 185.447876
##  [2,]  11.60801175 21.347420  30.089610
##  [3,]   0.07824755  3.323227 106.599051
##  [4,]  15.75374565 26.834280  24.948725
##  [5,]   9.15530094  3.247214  77.710675
##  [6,] 224.43179074 23.066347   4.313177
##  [7,]  15.28543337  7.866884   1.539904
##  [8,]  53.88994194 35.359990   7.051024
##  [9,]   7.21577099 15.631997  26.264064
## [10,]   1.52250306 25.688498  46.856173

ss.total <- sum(ss.especies)
ss.total

## [1] 1033.251

Para calcular os desvios quadráticos totais, nós subtraímos cada altura das árvores pela média geral e elevamos ao quadrado. Veja o gráfico:

vetor.cor <- rep(1:3, each=10) #vetor de cores

plot(x = c(1:30), y = arv3$alt, ylim=c(10,40),
     pch=(rep(c(15,16,17), each=10)),
     col=vetor.cor,
     ylab="Variável Resposta", xlab="Observações",
     main="Variação total")
    for(i in 1:30)
    {
    lines(c(i,i),c(arv3$alt[i],mean(arv3$alt)),col=vetor.cor[i])
    }
    abline(h=media.geral)

Agora vamos fazer a somatória dos desvios quadráticos dentro de cada grupo (ss.intra).

O gráfico para entender esse cálculo:

vetor.medias<-rep(medias.sps, each=10)

plot(c(1:30), arv3$alt, ylim=c(10,40),
     pch=(rep(c(15,16,17),each=10)),
     col=vetor.cor,
     main="Variação Intra Grupos",
     ylab="Variável Resposta", xlab="Observações")
    for(i in 1:30)
    {
    lines(c(i,i),c(vetor.medias[i],arv3$alt[i]),col=vetor.cor[i])
    }
    lines(c(1,10),c(medias.sps[1],medias.sps[1]),col=1)
    lines(c(11,20),c(medias.sps[2],medias.sps[2]),col=2)
    lines(c(21,30),c(medias.sps[3],medias.sps[3]),col=3)

Cálculo dos valores:

medias.sps

##      sp1      sp2      sp3 
## 20.44245 21.78334 27.62274

ss.sp1=sum((arv4$sp1-medias.sps["sp1"])^2)
ss.sp1

## [1] 262.2655

ss.sp2=sum((arv4$sp2-medias.sps["sp2"])^2)
ss.sp2

## [1] 157.0018

ss.sp3=sum((arv4$sp3-medias.sps["sp3"])^2)
ss.sp3

## [1] 322.4728

ss.intra=ss.sp1+ss.sp2+ss.sp3
ss.intra

## [1] 741.7401

A soma dos quadrados entre grupos:

plot(c(1:30), vetor.medias, ylim=c(10,40), 
  pch=(rep(c(15,16,17),each=10)),
  col=vetor.cor, 
  main="Variação Entre Grupos", 
  ylab="Variável Resposta", xlab="Observações")
 for(i in 1:30)
    {
    lines(c(i,i),c(vetor.medias[i],mean(vetor.medias)),col=vetor.cor[i])
    }
    abline(h=media.geral)
    points(c(1:30),arv3$alt, ylim=c(10,50), 
           pch=(rep(c(0,1,2),each=10)), col=vetor.cor, cex=0.5)

medias.sps

##      sp1      sp2      sp3 
## 20.44245 21.78334 27.62274

media.geral

## [1] 23.28284

ss.entre=10*sum((medias.sps-media.geral)^2)
ss.entre

## [1] 291.5112

#conferindo os cálculos
ss.intra+ss.entre

## [1] 1033.251

ss.total

## [1] 1033.251

Cálculo do F

# Desvios médios
ms.entre=ss.entre/2
ms.intra=ss.intra/27
ms.entre

## [1] 145.7556

ms.intra

## [1] 27.47186

#F - razão das variâncias
F.sps=ms.entre/ms.intra
F.sps

## [1] 5.305634

Distribuição F

Vamos ver na distribuição F qual a probabilidade de termos econtrado o valor F.sps ao acaso:

curve(expr=df(x, 2,27),
      main="Distribuição F de Fisher (df=2,27)", 
      xlab="Valor F",
      ylab="Densidade Probabilística (df)",
      xlim=c(2,12))
abline(v=F.sps, col="red")
abline(h=0, lty=2)

xf=seq(F.sps, 12, 0.01)
ydf=df(xf, 2, 27)
polygon(c(F.sps,xf),c(0,ydf),col="red")

text(x= 7,y=0.08,paste("pf(x) =",
     round(pf(F.sps,2,27,lower.tail=F),4)), 
     cex=0.8, col="red")

Cálculo do P

p.sps=pf(F.sps, 2, 27, lower.tail=FALSE)
p.sps

## [1] 0.01139248

Para mais informações sobre o F e as comparações com o teste t, veja esse roteiro.

A tabela final

Colocando os dados calculados em nossa tabela de anova

ANOVA no R

A função que constrói esta tabela de ANOVA e faz o teste estatístico é a aov. Vamos comparar:

anova.sps <- aov(alt~sp, data=arv3)
summary(anova.sps)

##             Df Sum Sq Mean Sq F value Pr(>F)  
## sp           2  291.5  145.76   5.306 0.0114 *
## Residuals   27  741.7   27.47                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Agora nós sabemos de ondem vem os números nesta tabela e como foi calculado o P, certo?

Onde estão as diferenças?

• veja os gráficos

• faça teste à posteriori - Tukey HDS

TukeyHSD(anova.sps)

##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = alt ~ sp, data = arv3)
## 
## $sp
##             diff        lwr       upr     p adj
## sp2-sp1 1.340893 -4.4708805  7.152667 0.8360276
## sp3-sp1 7.180300  1.3685259 12.992073 0.0132022
## sp3-sp2 5.839406  0.0276327 11.651180 0.0487468

#revendo gráfico dos dados
stripchart( arv3$alt~arv3$sp,
           vertical = TRUE, pch = 16, 
           col = c("black", "red", "green"))

#boxplot
boxplot(arv4)

lm.anova.sp <- lm(alt~sp, data=arv3)
summary.aov(lm.anova.sp) #mesmo que o summary do aov

##             Df Sum Sq Mean Sq F value Pr(>F)  
## sp           2  291.5  145.76   5.306 0.0114 *
## Residuals   27  741.7   27.47                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# mostra o resultado de forma um pouco diferente, consegue entender?
summary(lm.anova.sp) 

## 
## Call:
## lm(formula = alt ~ sp, data = arv3)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -12.1407  -3.0002  -0.0742   3.2719   9.2780 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   20.442      1.657  12.334 1.32e-12 ***
## spsp2          1.341      2.344   0.572  0.57202    
## spsp3          7.180      2.344   3.063  0.00492 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 5.241 on 27 degrees of freedom
## Multiple R-squared:  0.2821, Adjusted R-squared:  0.229 
## F-statistic: 5.306 on 2 and 27 DF,  p-value: 0.01139

• Y contínuo - X continuo

alt.nutr <- lm(alt ~ nutri, data = arv3)
summary(alt.nutr)

## 
## Call:
## lm(formula = alt ~ nutri, data = arv3)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -6.7176 -1.7171  0.0088  1.4124  9.2072 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 14.38711    1.05832  13.594 7.44e-14 ***
## nutri        0.58006    0.05971   9.714 1.82e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.906 on 28 degrees of freedom
## Multiple R-squared:  0.7712, Adjusted R-squared:  0.763 
## F-statistic: 94.37 on 1 and 28 DF,  p-value: 1.818e-10

• Y continuo

• X1 contínuo - X2 categórico

alt.sp.nutr <- lm(alt ~ nutri + sp, data = arv3)
summary(alt.sp.nutr)

## 
## Call:
## lm(formula = alt ~ nutri + sp, data = arv3)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.6682 -0.9554 -0.1585  1.3581  7.6134 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  13.3971     1.1034  12.142 3.23e-12 ***
## nutri         0.5256     0.0561   9.369 8.09e-10 ***
## spsp2         1.6421     1.1423   1.437  0.16251    
## spsp3         3.8326     1.1965   3.203  0.00357 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.553 on 26 degrees of freedom
## Multiple R-squared:  0.836,  Adjusted R-squared:  0.817 
## F-statistic: 44.16 on 3 and 26 DF,  p-value: 2.401e-10

alt.sp.nutr3 <- lm(alt ~ nutri*sp, data=arv3)
summary(alt.sp.nutr3)

## 
## Call:
## lm(formula = alt ~ nutri * sp, data = arv3)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.6724 -0.9444 -0.3313  1.0416  7.0702 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 13.40167    1.42727   9.390 1.66e-09 ***
## nutri        0.52527    0.08906   5.898 4.38e-06 ***
## spsp2        2.93717    1.97080   1.490    0.149    
## spsp3        0.43636    2.75666   0.158    0.876    
## nutri:spsp2 -0.10095    0.12423  -0.813    0.424    
## nutri:spsp3  0.17187    0.14350   1.198    0.243    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.474 on 24 degrees of freedom
## Multiple R-squared:  0.8578, Adjusted R-squared:  0.8282 
## F-statistic: 28.96 on 5 and 24 DF,  p-value: 2.003e-09

summary.aov(alt.sp.nutr3)

##             Df Sum Sq Mean Sq F value   Pr(>F)    
## nutri        1  796.8   796.8 130.179 3.53e-11 ***
## sp           2   66.9    33.5   5.467   0.0111 *  
## nutri:sp     2   22.6    11.3   1.846   0.1796    
## Residuals   24  146.9     6.1                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

• Outras distribuições de probabilidades além da Normal

• Variâncias não homogêneas

Mais comuns:

• Poisson – úteis para dados de contagem

• Binomial – úteis para dados com proporções, ou presença/ausência

1. Uma hipótese sobre a distribuição da variável resposta Yi.

2. Especificação do modelo (variáveis X).

3.Função de ligação entre a média Yi e as variáveis X. Função que lineariza a regressão

• Y - Número de atropelamentos de anuros em uma estrada

• X - distância a um Parque Nacional.

mod <- glm(TOT.N ~ D.PARK, data = atrop, family = poisson)
summary(mod)

## 
## Call:
## glm(formula = TOT.N ~ D.PARK, family = poisson, data = atrop)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -8.1100  -1.6950  -0.4708   1.4206   7.3337  
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  4.316e+00  4.322e-02   99.87   <2e-16 ***
## D.PARK      -1.059e-04  4.387e-06  -24.13   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 1071.4  on 51  degrees of freedom
## Residual deviance:  390.9  on 50  degrees of freedom
## AIC: 634.29
## 
## Number of Fisher Scoring iterations: 4

anova(mod, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: poisson, link: log
## 
## Response: TOT.N
## 
## Terms added sequentially (first to last)
## 
## 
##        Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
## NULL                      51     1071.4              
## D.PARK  1   680.55        50      390.9 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

##   (Intercept)        D.PARK 
##  4.3164849800 -0.0001058506

• Y - Quantidade de veados infectados dentre os capturados

• X - quantidade de vegetação aberta

nao.infectado <- veados$amostrado- veados$infectado
y <- cbind(veados$infectado, nao.infectado)
mod.veado <- glm(y ~  openland, data=veados, family="binomial")
summary(mod.veado)

## 
## Call:
## glm(formula = y ~ openland, family = "binomial", data = veados)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -3.8052  -1.1463   0.6707   1.3416   3.5483  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   2.4820     0.1736   14.29   <2e-16 ***
## openland     -6.0336     0.4432  -13.61   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 376.18  on 31  degrees of freedom
## Residual deviance: 121.30  on 30  degrees of freedom
## AIC: 204.65
## 
## Number of Fisher Scoring iterations: 5

anova(mod.veado, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: binomial, link: logit
## 
## Response: y
## 
## Terms added sequentially (first to last)
## 
## 
##          Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
## NULL                        31     376.18              
## openland  1   254.88        30     121.30 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# salvando os coeficientes
coefs <- coef(mod.veado)
coefs

## (Intercept)    openland 
##    2.482035   -6.033583

#lembre que estes coefs foram linearizados com a função de ligação
# para a binomial a link function é a logit
# então temos que transformar os y no inverso do logit, a função é essa:
invlogit = function(x) return(exp(x)/(1+exp(x)))

Teste de Bartlett:

bartlett.test(arv3$alt~arv3$sp)

## 
##  Bartlett test of homogeneity of variances
## 
## data:  arv3$alt by arv3$sp
## Bartlett's K-squared = 1.1109, df = 2, p-value = 0.5738

Teste de Levene:

leveneTest(arv3$alt~arv3$sp)

## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  2   0.512  0.605
##       27

A premissa de normalidade na verdade recai sobre os resíduos do modelo.

resid(anova.sps)

##           1           6          11          16          21          26 
##   0.8395548  -0.5666556   3.1201253  -1.1287012   5.8661704 -12.1406501 
##          31          36          41          46          51          56 
##   6.7500566  -4.5005793   0.1541789   1.6065001  -2.6382600  -3.1208224 
##          61          66          71          76          81          86 
##   3.3224765  -3.6806771  -0.3024985   6.3022481  -1.3052922   7.4459310 
##          91          96         101         106         111         116 
##  -2.4542276  -3.5688778   9.2780228   1.1454976   5.9847796  -9.3347719 
##         121         126         131         136         141         146 
##   4.4754638  -2.2630830  -5.5808309  -6.9952785   0.7849457   2.5052549

#residuals(anova.sps) # o mesmo que a função de cima